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Quantum J1-J2 antiferromagnet on the stacked square lattice: Influence of the 
interlayer couphng on the ground-state magnetic ordering 

D. Schmalfufi,^ R. Darradi,^ J. Richter,^ J. Schulenburg,^ and D. Ihle^ 

^ Institut fiir Theoretische Physik, 

Umversitdt Magdeburg, 39016 Magdeburg, Germany 
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Using the coupled-cluster method (CCM) and the rotation-invariant Green's function method 
(RGM), we study the influence of the interlayer coupling J± on the magnetic ordering in the ground 
state of the spin-1/2 J1-J2 frustrated Heisenberg antiferromagnet (J1-J2 model) on the stacked 
square lattice. In agreement with known results for the J1-J2 model on the strictly two-dimensional 
square lattice {J± =0) we find that the phases with magnetic long-range order at small J2 < Jci and 
large J2 > Jc2 are separated by a magnetically disordered (quantum paramagnetic) ground-state 
phase. Increasing the interlayer coupling J± > the parameter region of this phase decreases, and, 
finally, the quantum paramagnetic phase disappears for quite small J± ~ 0.2 — 0.3Ji. 



The properties of the frustrated spin-1/2 Heisen- 
berg antiferromagnet (HAFM) with nearest-neighbor Ji 
and competing next-nearest-neighbor J2 couphng ( J1-J2 
model) on the square lattice have attracted a great deal of 
interest during the last fifteen years (see, e.g., Refs. 1-12 
and references therein) . The recent synthesis of layered 
magnetic materials^'^'^'* which can be described by the 
J1-J2 model has stimulated a renewed interest in this 
model. It is well-accepted that the model exhibits two 
magnetically long-range ordered phases at small and at 
large J2 separated by an intermediate quantum param- 
agnetic phase without magnetic long-range order (LRO) 
in the parameter region Jc^ < J2 < Jc2j where J^ « 0.4 
and Jc2 ~ 0.6. The ground state (GS) at low J2 < J^ 
exhibits semi-classical Neel magnetic LRO with the mag- 
netic wave vector Qq = (tTjTt). The GS at large J2 > Jc2 
shows so-called collinear magnetic LRO with the mag- 
netic wave vectors Qi = (7r,0) or Q2 = (0, tt). These 
two collinear states are characterized by a parallel spin 
orientation of nearest neighbors in vertical (horizontal) 
direction and an antiparallel spin orientation of nearest 
neighbors in horizontal (vertical) direction. The prop- 
erties of the intermediate quantum paramagnetic phase 
are still under discussion, however, a valence-bond crys- 
tal phase seems to be most favorable. ^"'^'^'^ 

The properties of quantum magnets strongly depend 
on the dimensionality.^^ Though the tendency to order is 
more pronounced in three-dimensional (3d) systems than 
in low-dimensional ones, a magnetically disordered phase 
can also be observed in frustrated 3d systems such as 
the HAFM on the pyrochlore lattice^^ or on the stacked 
kagome lattice. ^^ On the other hand, recently it has been 
found that the 3d J1-J2 model on the body-centered cu- 
bic lattice does not have an intermediate quantum para- 
magnetic phase. ^^'^^ Moreover, in experimental realiza- 
tions of the J1-J2 model the magnetic couplings are ex- 
pected to be not strictly 2d, but a finite interlayer cou- 
pling Jx is present. For example, recently Rosner et 



al.^^ have found Jj_/Ji ^ 0.07 for Li2VOSi04, a ma- 
terial which can be described by a square lattice J1-J2 
model with large J2-^^'^'^ 

This motivates us to consider an extension of the J1-J2 
model, namely the J1-J2 spin-1/2 HAFM on the stacked 
square lattice described by the Hamiltonian 

H = 2_^ I "^1 /_^ ^i,n ■ Sj,n + J2 2_^ ^i,n ' ^j.n I 



-J-T.'- 



^i.n-\-l ; 



(1) 



where n labels the layers and J^ > is the interlayer 
coupling. The expression in brackets represents the Ji- J2 
model of the layer n with intralayer couplings Ji = 1 and 
J2 > 0. The main problem we would like to study is the 
influence of J± on the existence of the intermediate quan- 
tum paramagnetic GS phase. Note that the exact diago- 
nalization widely used for the study of the strictly 2d Ji- 
J2 model, see, e.g., Refs. 2-4, is not appropriate for the 
3d problem under consideration. Therefore, we use the 
coupled-cluster method (CCM)^'^°"^^ and the rotation- 
invariant Green's function method (RGM). '^'^'^■^s. 28-30 
Both methods have been successfully applied to quan- 
tum spin systems in arbitrary dimension and are able to 
deal with frustration. 

Let us briefly illustrate some basic features of the 
CCM. For more details the reader is referred to Refs. 
6 and 20-24. The starting point for the CCM calcula- 
tion is the choice of a reference state |<i>). For |(f>) of 
the considered spin system we choose the two-sublattice 
Neel state for small J2 but a collinear state for large J2. 
To treat each site equivalently we perform a rotation of 
the local axis of the spins such that all spins in the ref- 
erence state align along the negative z axis, i.e., in the 
rotated coordinate frame we have |$) = | J.)| i)| i) • . . • 
Note that in this new frame the Hamiltonian is modified. 



and |*&)=|i)|J,)|i) ... is not an eigenstate of this modified 
Hamiltonian, see, e.g. Refs. 20,21,23. For the kct GS |*) 
with H\'^) = E\'^) an exponential ansatz |$) = 



is used, where the correlation operator S is given by 
S ~ ^i-io'^i^i ■ "^^^ ^i represent a set of muhi-spin 
creation operators Cj = sj' , s^s 



'^ ""j 'k 



The application of all the Cf on |$) creates a complete 
set of states, which may contribute to |\E'). The corre- 
lation operator S contains the coefficients Sj which are 
determined by requiring that the expectation value of H 
is a minimum. The order parameter M is given by the 
expectation value of s^ . 

For the considered quantum many-body model it is 
necessary to use approximations in order to truncate 
the expansion of 5*. We use the well elaborated LSUBn 
scheme^^'^^'^"^ in which in the correlation operator S all 
multi-spin correlations over all distinct locales on the lat- 
tice defined by n or fewer contiguous sites are taken into 
account. For example, within the LSUB4 scheme one in- 
cludes multi-spin creation operators of one, two, three or 
four spins distributed on arbitrary clusters of four con- 
tiguous lattice sites. The number of these fundamental 
configurations can be reduced exploiting lattice symme- 
try and conservation laws. In the CCM-LSUB8 approx- 
imation we have finally 25953 (43070) fundamental con- 
figurations for the Neel (coUinear) reference state. To 
solve the set of the corresponding ket equations we use 
parallel computing. ^^ 

Since the LSUBn approximation becomes exact for 
n —> oo, it is useful to extrapolate the 'raw' LSUBn data 
to n ^ oo. An appropriate extrapolation rule for the or- 
der parameter of systems showing a GS order-disorder 
transition is the 'leading power-law' extrapolation^"^ 
M{n) = Co + ci{l/nY^, where the results of the 
LSUB4,6,8 approximations are used for the extrapola- 
tion. For the GS energy per spin e{n) — ao + ai(l/ri^) + 
a2{l/n'^) is a reasonable extrapolation ansatz. ^^ 

Next we give a brief illustration of the spin-rotation- 
invariant Green's function method. ^^'^^ More details can 
be found in Refs. 5,17,28 and 30. Considering the 
equations of motion for the commutator Green's func- 
tion (\Sq ; Sq )) and supposing spin-rotation invariance, 
i.e. (4J = 0, we get uj^ {{s+; s')) ^ = ([js+,Sq]_) + 

((— Sq ; Sq )) . To treat the operator Sq containing prod- 
ucts of three spin operators along nearest-neighbor se- 
quences, a decoupling procedure in the spirit of Ref. 26 
is performed. For example, the operator product s^sj^s'^ 

is replaced by tja^b (•s^^b) ^c + ^A,c (s^^c) *B' where 
A, B, C represent spin sites. The introduction of vertex 
parameters ry-y^^ is aimed to improve the approximation 
and to fulfill fundamental constraints like the sum rule. 
By analogy with Refs. 5 and 28 we use four differ- 
ent vertex parameters, namely ?7i|| related to the corre- 
lator ci^o.O: 'yix related to co.o,i, ?72 commonly related 
to C2,o,o, C2,i,o, C2,2,o, Ci,04, Ci^^i, co,o,2, and 773 re- 
lated to ci^iQ. The correlators are defined as Ck,i,m = 
cr = (sq Sr) — 2 (sqSr) /3 with the lattice vector 



R = fcai -I- la2 + maa and have to be determined self- 
consistently. Performing the approximations mentioned 



1$) above we obtain ((s+;Sq)) — 



mq/{Lo'^ — cjq), where 



for 



and cijq explicit equations can be given. The 



equation for u'i contains the four vertex parameters and 
the nine correlators mentioned above. The correlators 
can be expressed by the Green's function using the spec- 
tral theorem. To determine the four vertex parameters 
we use the sum rule cq^o.o = 1/2 and require that the 
static susceptibility x+^ = — lim^^^o ((sq ; Sq)) has 
to be isotropic in the limit q -^ q i7,28,30 rpj^g j-gj^jain- 
ing two equations arc obtained as follows: First we use 
the relation 773 = (7726^"'^ + J2ili\\){^ + ^2)^^ which was 
successfully applied in Ref. 5 to the 2d J1-J2 model. 
This relation interpolates between the two limiting cases 
J2 ^ and J2 ^ 00 and takes care of the relation 
limj^^o ci,o,o — fin:ij2-»oo ci,i,o. Finally we use, following 
Ref. 5, an approximative expression for the GS energy per 
spin e™''" = 3 JiCi^o, 0+3 J2C1, 1,0+3 Jj_ Co, 0,1/2 as an addi- 
tional input. For the stacked HAFM considered we make 
the ansatz eo"^"*(J2, </j_) = fi{Ji.) + /2(>/2) (note that 
Ji = 1). To fix /2 we use the exact diagonalization result 
for the GS energy of the finite 2d Ji- J2 model ( J^ = 0) 



of iV = 32 spins, i.e. we set ,f2{J2) 



^N=32 



{J2,J±^0). 



To fix /i we use the GS energy of the unfrustrated 
stacked square lattice Cq^ {J2 — 0,J±) calculated by 
linear spin- wave theory and set fi{J±) = e^^ {J2 = 
0, J±) — eQ^^'^{J2 = 0,J± = 0) this way taking into ac- 
count the effect of the intcrlayer coupling and a finite-size 
correction. 

To discuss GS magnetic order-disorder transitions we 
consider the magnetic order parameter. In the RGM 
scheme^^"^^'^° the correlation function (sqSr) at T = 
is given by 



(sqSr) 



3 

27V 



"^ 2wa 



Qi 



Qi 



-iQjR 



(2) 



The second term (condensation part) describes LRO, 
where the sum runs over different noncquivalcnt mag- 
netic wave vectors Qj taking into account the possi- 
bility to have degenerate GSs. For model (1) we have 
Qo = (7r,7r,7r) for the Neel phase and Qi — (7r,0,7r) or 
Q2 = (0,7r,7r) for the coUinear phase. Magnetic LRO 
is accompanied by a diverging static susceptibility Xq~ 
at q = Q giving an additional equation for Cq. Note 
that for the coUinear phase both condensation terms are 
equal, i.e. Cq^ ~ Cq^. The order parameter M can be 
calculated by M^ — 3 |Cq| /2. That way, the order pa- 
rameter is linked to the long-range behavior of the corre- 
lation functions because M is nonzero if limiRi^oo (sqSr) 
remains finite. 

As in the 2d case the GS of the stacked model is char- 
acterized by two magnetically long-range ordered phases, 
namely a Neel phase for small J2 and a coUinear phase for 
large J2. For not too large J^ both magnetic phases are 
separated by a magnetically disordered quantum param- 
agnetic phase, where the phase transition points arc func- 
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FIG. 1: Magnetic order parameter M versus J2 for various 
strengths of the interlayer coupling Jx . (a): RGM, (b): CCM 



tions of J ^ . To determine these transition points we cal- 
culate the order parameters for various Jj^ to find those 
values J2 — OLNeeiiJi.) and J2 — acouiJ±) where the or- 
der parameters vanish. In Fig. 1 we present some typical 
curves showing the order parameters versus J2 for some 
values of J± . Both approaches lead to qualitatively com- 
parable results. The magnetic order parameters of both 
magnetically long-range ordered phases vanish continu- 
ously as it is typical for second-order transitions. Note, 
however, that there are arguments'^'^ that the transition 
from the collinear phase to the quantum paramagnetic 
phase should be of first order. The order parameters 
are monotonously increasing with J± , and the transition 
points aNeei and acoii niove together. In Fig. 2 we present 
these transition points in dependence on J± . Close to the 
strictly 2d case, i.e. for small J_l ^ 1, the infiuence of the 
interlayer coupling is largest. For a characteristic value 
of Jl w 0.31 (0.19) for the RGM (CCM) approach the 
transition points ajsieei and Ucoii meet each other. 

For larger J^ exceeding J\_ we have a direct first-order 
transition between both types of magnetic LRO as it is 
also observed in the classical model and in the 3d quan- 
tum J1-J2 model on the body-centered cubic lattice. ^^'^^ 



FIG. 2: Ground-state phase diagram, (a): RGM, (b): CCM. 
The solid lines show those values of J2 where the order param- 
eters vanish. The dashed line in (b) represents those values of 
J2 where the two energies calculated for the Neel and collinear 
reference state become equal. 



However, the description of this first-order transition is 
not possible within the RGM approach. The reason is 
that the approximative expression for the GS energy per 
spin e™^" used as an input is a smooth function of J2, 
whereas a first-order GS transition is characterized by a 
kink in gq. As a consequence we find that there is no 
solution of the system of coupled RGM equations for pa- 
rameter values being close to a first-order transition, i.e. 



for J2 ~ 0.5 and J± > J*^. The order parameter curve 
for J_L = 0.4 depicted in Fig. 1(a) indeed shows a small 
region slightly below J2 = 0.5, where no solution exists. 
In contrast to the RGM the CCM approach starts with 
two different reference states (Neel and collinear) related 
to the two types of magnetic LRO. Though we start our 
CCM calculation with a reference state corresponding to 
semiclassical order, one can compute the GS energy also 
in parameter regions where semiclassical magnetic LRO 
is destroyed, and it is known^'^^"^^ that the CCM yields 
precise results for the GS energy beyond the transition 
from the semiclassical magnetic phase to the quantum 
paramagnetic phase. The necessary condition for the 
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FIG. 3: CCM results for the energy per spin e for both refer- 
ence states (inset) and the order parameter M for Jx = 0.2. 
Both quantities are obtained by extrapolation of the 'raw' 
LSUBn results to the limit n — > oo as explained in the 
text. The energies calculated with Neel and coUinear refer- 
ence states become equal at J2 « 0.58 indicating a first-order 
transition. For the order parameter M we take that value 
calculated with the reference state of lower CCM energy. 



convergence of the CCM equations is a sufficient overlap 
between the reference state and the true GS. Hence we 
can add to the above discussion of the order parameters 
a comparison of the energies. Provided that the CCM 
equations converge for the Neel and the coUincar refer- 



ence state far enough beyond those points where the or- 
der parameters vanish we can determine the point where 
both energies become equal. For the considered LSUBn 
approximations this happens for J± ^ 0.1. In the inset 
of Fig. 3 wc show the energies versus Ji for Jj^ = 0.2 
caclulated by extrapolation. The corresponding points 
J2 = ctJ,Qj;(Jj^) where both energies meet are shown in 
Fig. 2 as dashed line. 

We obtain that both transition points UcoU and a^^;; 
are close to each other and show a similar dependence 
on J_L. Secondly, we find that at least for J^ > 0.1 the 
energy obtained with the Neel reference state is lower 
than that obtained with the coUinear reference state 
even for J2 values where the Neel order parameter is al- 
ready zero but the coUinear order parameter is still finite. 
Thus, this energetic consideration leads to the following 
sequence of zero-temperature transitions: Second-order 
transition from Neel LRO to a quantum paramagnetic 
phase at J2 = otNeei and then a first-order transition 
from the quantum paramagnetic phase to coUinear LRO 



at J2 



a. 



> Oicoii > ctNeei- This behavior is illus- 



trated in Fig. 3, where the order parameter M is shown 
versus J2 for fixed J± — 0.2. For a certain value of 
J± « 0.23 both transition points aNeei and a'^^n become 
equal, and one has a direct first-order transition between 
the two semiclassically long-range ordered phases. 
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